Section: New Results

Supercomputing for Helmholtz problems

The numerical solution of Helmholtz problems set in highly heterogeneous media is a tricky task. Classical high order discretizations fail to handle such propagation media, because they are not able to capture any of the scales of the velocity parameter. Indeed, they are build upon coarse meshes and therefore, if the velocity parameter is taken to be constant in each cell (through averaging, or local homogenization strategy), scale information is (at least partially) lost. We propose to overcome this difficulty by introducing a multiscale medium approximation strategy. The velocity parameter is not assumed to be constant on each cell, but on a submesh of each cell. If the submeshes are designed properly, the medium approximation method is equivalent to a quadrature formula, adapted to the medium. In particular, we show that this methodology has roughly the same computational cost as the classical finite element method. This new solution methodology has been presented in a paper under revision. We have performed a mathematical analysis of the multiscale medium approximation techniques to higher order discretization. First, we show that the heterogeneous Helmholtz problem is well-posed and derive stability estimates with respect to the right hand side, and with respect to variations of the velocity parameter, justifying the use of medium approximation. Those results are obtained assuming the velocity parameter is monotonous and that the propagation medium is closed by first order absorbing boundary conditions. However, these hypothesis are not mandatory to discretize the problem. Second, we turn to the analysis of finite element schemes with subcell variations of the velocity. In particular, we show that even if the solution can be rough inside each cell because of velocity jumps, we are able to extend the asymptotic error estimates obtained in  [93] to heterogeneous media with non-matching mesh in case of elements of order $1\le p\le 3$. Third, we investigate numerically the stability of the scheme when the frequency is increasing to figure out optimal meshing conditions. We show that in simple media, the optimal homogeneous pre-asymptotic error estimates are still valid. However, in more complex cases, it looks like this condition is not sufficient anymore. Apart from showing that the homogeneous results are not always applicable to the heterogeneous Helmholtz equation, we are not able to give a clear answer to the question. Finally, we are able to conclude that high order methods are actually interesting: in our examples, $p=4$ discretizations always yield a smaller linear system than lower order discretizations for the same precision.

For the numerical analysis of the scheme, we have shown that the HDG method could be rewritten as an upwind fluxes DG method and one of our perspectives is to use this equivalence in order to perform a dispersion analysis following the work of Ainsworth, Monk and Muniz  [64] .

We have shown that HDG could be used for 2D simulation on geophysical benchmark, and we will now implement the method in a Reverse Time Migration software, the ultimate goal being to couple HDG method with a full waveform inversion solver. In order to tackle more realistic test cases in 3D, it will be mandatory to improve the linear solver and we are now considering the use of an hybrid solver such as Maphys developed by the Inria team-project HIEPACS.

The results of this work have been presented at the “SIAM Conference on Geosciences” [48] and at the “Oil and Gas HPC Workshop” [49] .